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Abstract 

The small change in the spectral slope of the overall intensity of cosmic rays near 1 
PeV may be associated with the endpoint energy of supernova shock acceleration. A 
crucial test of this connection and other ideas of the origin of the spectral 'knee' is 
the reliable determination of the variation of elemental composition in this region. 
Recent measurements at the DICE/CASA/MIA air shower installation in Dugway, 
Utah, USA have provided several independent air shower parameters for each event. 
The derivation of elemental composition from a combination of Cherenkov size, 
depth of shower maximum in the atmosphere, muon size and electron size at ground 
level and the reliability of these results are discussed. There is no evidence from these 
data for a large change in the mean mass of cosmic rays across the 'knee'. These 
measurements show the cosmic rays are composed of ^ 70% protons and a particles 
near total energies of lOPeV. 



1 Introduction 



The slight steepening in the power law spectrum of the cosmic ray flux near 
energies of 10^^'^ eV has provoked a wide range of interpretations since its 
discovery in air shower size spectra nearly 40 years ago[l]. This spectral 'knee' 
remains an enigma, which together with the more recently discovered flat- 
tening ( or 'ankle') near 10^^'^ eV provide the only known features in the 
arriving cosmic ray spectrum at high energy. The coincidence of the 'knee' 
with the maximum energies expected from diffusive shock acceleration in su- 
pernova remnants, has served to cement its value as a key region for further 
investigation in cosmic-ray research. The features which might better connect 
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the change in overall spectral slope with expectations from shock acceleration 
lie in the details of the individual elemental spectra hidden within the over- 
all cosmic-ray flux. The accurate determination of the elemental composition 
across the knee region, from ~ 10^^ to lO^^eV, is a crucial measurement yet 
to be made. 

Recent observations of supernovae using TeV 7 rays have not managed to de- 
tect [2] the signal expected from high-energy cosmic-ray collisions with local 
material if supernova remnants (SNR) are the sources of high-energy cosmic 
rays [3] . Although EGRET has reported Ge V 7 rays from some SNRs [4] , it is 
still unclear if the signal is produced by locally intense cosmic rays [5]. Non- 
thermal x-ray emissions have been detected from two SNRs [6,7] which have 
been interpreted as synchrotron emission from local electrons at energies up 
to ~100 TeV [8], but the source spectral slope inferred by these results seems 
inconsistent with the electron source spectrum derived from direct measure- 
ments. It remains the case that the most convincing argument for supernovae 
explosions producing the bulk of the cosmic rays arises from the staggering 
power required to maintain the population given the measured lifetime of par- 
ticles within our Galaxy. A typical value for the power needed is ~ 10^°erg/s 
(see e.g. Schhckeiser[9]), based on a cosmic-ray residence time of ~ 10^ years, a 
time-scale directly inferred from measurements of the surviving amount of the 
radioactive isotopes of Be^^'flO] and AP^[11] which are created by spallation 
reactions of cosmic rays on interstellar material. 



2 Air Shower Detectors 

While direct detection above the atmosphere is clearly the most desirable 
method for composition determination, the fluxes of particles in this energy 
range are so small that a detector large enough to collect a good statistical 
sample up to lOPeV (lPeV=10^^eV) has not yet been built. Indirect detec- 
tion of cosmic rays, through measurements of air showers produced in the 
atmosphere, can easily supply sufficient collecting power but these methods 
generally have poor resolution for the properties of the originating cosmic ray 
nucleus whose identity becomes clouded by the huge number of interactions 
in the air shower. Further, the reconstruction of the mass of the incoming nu- 
cleus from measurements of the shower distributions at ground level is subject 
to systematic errors introduced by reliance on shower modeling by numerical 
simulations. 

To address these issues we have constructed and operated a new type of exper- 
iment designed specifically for cosmic-ray measurements near the knee. This 
Dual Imaging Cherenkov Experiment (DICE) is a ground based air shower 
detector which is designed to have as httle reliance as possible on the details 
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of the air shower simulations and to have the capabihty of comparison with 
existing direct measurements at O.lPeV to provide an assessment of the overall 
systematic error. DICE operates at the site of the Chicago Air Shower Array 
(CASA) and the Michigan Array (MIA) which are designed as an observatory 
for PeV 7 rays. They can also be used to examine hadronic showers since they 
give measurements of electron size and muon size at ground level. 

2.1 DICE description 

Since the method of imaging hadronic showers in Cherenkov light is a relatively 
recent development, we provide some more detailed description of the DICE 
detectors and operation. The two DICE telescopes are located at the CASA- 
MIA site in Dugway, Utah (described in Borione et al. (1994)[12]). They each 
consist of a 2m diameter f/1.16 spherical mirror with a focal plane detector 
of 256 close packed 40mm hexagonal photomultipliers (PMTs) which provide 
~ 1° pixels in an overall field of view 16° x 13.5° centered about the vertical. 
The telescopes are on fixed mounts separated by 100m, (see Boothby et al. 
(1995, 1997) [13,14]). 

Cosmic-ray events within the field of view produce a focal plane image at the 
PMT cluster which corresponds to the direction and intensity of Cherenkov 
light coming from the air shower. When the direction of the air shower core and 
the distance of the shower axis from the telescopes are known, simple geometry 
can be used to reconstruct the amount of light received from each altitude of 
the shower. The amount of Cherenkov light produced is strongly correlated 
with the number of electrons in the shower and is used to estimate the electron 
size as a function of depth in the atmosphere from which the location of shower 
maximum can be determined. This procedure is essentially geometrical and is 
independent of numerical simulations except for calculations which determine 
the angular distribution of Cherenkov light around the shower axis. 

2.2 CAS A/MIA description 

The CASA/MIA installation and operation is described in detail in Borione 
et al. 1994[12]. In this present work, information from CASA-MIA is used to 
establish the event geometry and to directly measure the electron size and 
muon size at ground level. Since the event acceptance is constrained by the 
DICE aperture, only events within ~10° of the zenith are analyzed, making 
the atmospheric depth of the size measurements to be essentially the same for 
all events at ~860g/cm^. 
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Table 1 

Parameters measured for each event. 



Description 


Symbol used here 


Source Instrument 


Shower direction and core location 




CASA 


Chcrcnkov size 


Ch 


DICE 


Shower maximum 




DICE 


Electron size at ground 




CASA 


Muon size at ground 




MIA 



3 Measured Parameters 



For each air shower collected a simple time coincidence is used to identify the 
same event in both DICE clusters and CASA-MIA. Further requirements on 
the correlation of the DICE images with each other and with the CASA-MIA 

event geometry arc Tiscd to reduce the overall probability of event mismatches 
between the detectors to ~10~^. The basic independent parameters measured 
for each shower are shown in Table 1 

The parameters for each shower are derived from these measured values. The 
accuracy of the shower core location derived by CASA is l-3m depending on 
the overall shower size. The measurement of the shower arrival direction is 
accurate to ~ 0.4" for larger showers with some degradation for lower en- 
ergy events. These measurements and their errors are discussed in detail in 
Borione et al. 1994[12]. The electron size from CASA is derived from a fit 
to the measured particle intensities to signals produced in an array of plas- 
tic scintillation counters using the NKG formulism[15,16]. The muon size is a 
fit to signals produced in buried scintillators using the function discussed by 
Greisen[15]. Further details of these fits are discussed in McKay[17]. 

3. 1 Derived Shower Parameters 

Some of the shower parameters used in this analysis correspond directly to 
measured values. The electron and muon sizes at ground level {Ne and Nf^) 
are those derived by CASA/MIA. The Cherenkov size at the two DICE de- 
tectors is obtained by summing the total amount of light detected at each 
photomultiplier tube (PMT) cluster. The location of shower maximum in the 
atmosphere (Xmax) is determined by fitting the shape of the shower image 
in each of the DICE PMT clusters. The Cherenkov light from each height in 
the shower falls into a specific angular bin in the focal plane when the event 
geometry is known. By applying corrections for solid angle, the variation of 
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light emission intensity with respect to the shower axis, and corrections for 
atmospheric absorption, a profile of electron size versus atmospheric depth is 
developed from the Cherenkov image. This is fit with a simple two parameter 
shower shape to find the location of shower maximum in each DICE telescope. 
An estimate for {X^ax) is found for both DICE clusters, a crucial difference 
between DICE and detectors which fit the lateral distribution of Cherenkov 
light is that this estimate of (Xmax) is based on simple geometry. It does 
not rely on the details of numerical simulations or have stringent linearity 
requirements of the detectors. 

Previous work with DICE estimated the shower energy by a simple transla- 
tion from the total amount of Cherenkov light in the image and the geometry 
of the shower [18]. In this present work a more accurate estimate of energy is 
derived from a combination of the amount of Cherenkov light and the X^ax 
determination produced by each DICE telescope. This is desirable since the 
lateral distribution and intensity of Cherenkov light at a given total energy de- 
pends both on the primary particle mass, hence mean Xmax, and the distance 
of the measurement from the shower core. A fit for the total shower energy 
and primary particle mass is made to the geometry, Cherenkov size {Ch) and 
Xmax location in the two DICE telescopes. The form of the Cherenkov size 
function used in these fits is derived from the results of simulations using the 
program CORSIKA 4.50 with the VENUS interaction model([19]). The de- 
rived Xmax fitting function has a constant shower elongation rate parameter 
and assumes simple superposition for providing primary particle mass depen- 
dence. The exact form of both of these fitting functions are given in Appendix 
A as equations A.l and A. 2. 

Events collected by DICE in coincidence with CASA/MIA over a period from 
mid 1994 to early 1996 are subjected to the following selection cuts before fit- 
ting with equations A.l and A. 2. (1) The core of the shower lies at a distance 
100m < r < 225r/), from both DICE telescopes. (2) The fits of the longitudinal 
development in both DICE telescopes have reduced < 3. (3) The Xmax 
from each telescope agrees within IbOg/cm^. (4) The arrival direction of the 
shower is within 6° of the vertical. The fits to the measured Ch and Xmax 
from the two sites are made using the minimum method. With these cuts 
the effective collection geometry is determined by the instrument Monte Carlo 
to be ~3300 m^sr, making the overall collecting power ^125,000 m^sr days. 
Using the detailed calculated aperture and efficiency corrections derived from 
an instrument Monte Carlo discussed below the cosmic ray energy spectrum 
can be constructed from the fit energies. This is shown as the filled circles in 
Figure 1 together with previous determinations of the energy spectra in this 
region. The fiuxes have been multiplied by i?^-^^ to emphasize the change in 
slope across the knee region. These agree with the previously reported energy 
spectrum for DICE [18]. We believe this new analysis provides a better deter- 
mination of fluxes because of the primary mass dependent fltting function used 
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Fig. 1. Measurements of energy spectra near the 'knee' from this work (filled circles) 
compared with recent measurements from Tibct[23], CASA-MIA [24], HEGRA[20], 
Tunka[21], and Akeno[25]. 



for each event. These new DICE data are given in tabular form in Appendix 
B. The 'knee' here is at an energy around 3PeV. 

The mean measured values of Xmax as a function of energy is shown in Figure 
2. Here these values are compared with the previous DICE results[18] and 
those of the HEGRA experiment determination [20] which uses the distribu- 
tion of Cherenkov light at ground level to estimate X^ax- The star shows the 
expectation for the mean X^ax based on a the direct composition measured 
near 100TeV[22] using equation A. 2. The various measurements show general 
agreement, although there is a tendency for the DICE results to show a lighter 
composition at the highest energy. The present data are in statistical agree- 
ment with the previously published results of DICE, but they show somewhat 
less departure from constant composition at high energies. The dashed hues 
are for pure proton or pure iron composition given by equation A. 2. These 
results are tabulated in Appendix B. 
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Fig. 2. Measurements of mean X, 
3.2 Robustness of Shower Parameters 
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The determination of incident particle mass from air shower measurements 
have generally poor resolution. Consequently features such as elemental peaks, 
which can commonly be used in direct measurements of cosmic rays to test 
detector calibration, are absent. This makes systematic sources of error the 
largest obstacle which must be overcome in order to produce a reliable mea- 
surement. Any internal consistency checks which can be used on the data are 
valuable as a way to keep these problems to a minimum. The combination 
of DICE and CASA/MIA provides a unique combination of shower parame- 
ters which over-determine the characteristics of each shower. The comparison 
of parameters from different detectors can be used to test the rehabihty of 
results. This process is already possible within DICE since the separate mea- 
surements of Xmax and Ch made with each cluster can be compared to test 
consistency. A more rigorous test is a comparison of the information derived 
from DICE with that from CASA/MIA. 

The general profile of air showers produced by hadrons has been studied by 
several authors. These showers can be thought of as a succession of electro- 
magnetic showers fed by a core of interacting hadrons. The fluctuations in 
shower maximum are produced both by the statistical nature of the shower 
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development and the fluctuations in the first interaction point in the atmo- 
sphere. At atmospheric depths deeper than X^ax, the shower size dechnes in 
a manner which is dominated by the atmospheric hadronic interaction length. 
If the size is measured at a depth relative to the location of X^nax the fluctua- 
tions produced by variations in the initial interaction point are removed. We 
can perform a simple test on the data discussed here: Does the electron size 
at ground level measured by CASA fit with the shower energy and location of 
^max produced by DICE? By using the average shower development function 
quoted by Gaisser[26] (given in Appendix A), we can compare the expected 
ground electron size with the CASA measurement. 

Since showers at a given energy have values of Xmax which fiuctuate for the 
reasons outlined above, we can sample a range of potential electron sizes versus 
^max for a specific primary energy range. This procedure can be used to 
verify the consistency of measurements and to test for systematic errors in 
the determination of Xj^ax, E, and Ne as a function of energy. If at a ground 
level depth of Xq we compare Ne{Xo)/E with Xmax for each event we can 
test on average if these agree with equation A. 3 over a range of parameters. 
Figure 3 shows how this test works for Xj^ax and E determined from DICE 
and A^e from CASA. The four panels correspond to some of the energy bins 
used for Figure 2. Each panel shows the average loglO(Ae/-E') plotted versus 
shower Xmax- The function derived from the average air shower development 
(from equation A. 3) is shown as a solid line. The good agreement between the 
measurements and the data for the lower energy data provide confidence that 
the location of Xmax determined by DICE and the electron size at ground 
level are measured correctly. The version of the Gaisser function shown as the 
dashed curve in the lower right hand panel of Figure 3 shows the effect of 
scaling violation as discussed in Appendix A. 3. 

It is important to realize that not just the shower shape but the normalization 
of the shower size versus Xmax is given by equation A. 3, there is no arbitrary 
vertical normalization of the curves shown in Figure 3. These results therefore 
also provide confidence that the measurement of energy from DICE is reason- 
able, an error in energy would produce a systematic shift of the data in the 
vertical direction. For the higher energies there appears to be a small system- 
atic shift, of the order of 10%, of the data to a region below the curve. This 
could be produced by a various effects; (i) the energy scale is overestimated at 
the higher energies, (ii) the electron size is underestimated for large electron 
sizes, (iii) the location of Xmax is overestimated for high energies. This last 
possibility seems unlikely since the method used for the measurement of Xmax 
is the same at both low and high energies and no similar systematic shift ap- 
pears across all energies. For example the region for Xmaj; near 700g/cm^ is 
sampled both by the IPeV and lOPeV energy, at the lower energy the mea- 
surements agree with the expected curve, at higher energies they do not. A 
systematic error in the determination of Xmax for deep showers should pro- 



8 




xmax (g/cm ) xmax (g/cm ) 

Fig. 3. Comparisons of expected electron size at ground with a shower model from 
Gaisser [26]. The open symbols represent raw electron sizes, the filled symbols are 
electron sizes corrected for saturation in a manner discussed in the text. 

duce a similar shift at both energies. The lack of a systematic shift across all 
energies also shows that errors in X^ax produced by atmospheric scattering 
of Cherenkov light is minimal, since this would also only depend on the depth 
of maximum shower size. The first possibility, (i) mentioned above, also seems 
unlikely since there would be a large effect on the energy spectrum shown in 
Figure 1. A shift in the energy scale measure of the required magnitude (~ 
20%) would produce a significant increase in the apparent fluxes at the higher 
energies (~ 40%), which does not seem to be present in the data. We are left 
with the possibility that the electron size somehow 'saturates' slightly at the 
higher energies if we believe equation A. 3 is independent of energy. This is a 
possibility since CASA, primarily built for observations near lOOTeV, was not 
designed to operate for air shower core density measurements at lOPeV. The 
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electron sizes at high energy are predominantly derived from counters on the 
perimeter of the core which have not saturated. We have developed a simple 
correction to the electron size for large values which apparently corrects the 
discrepancy shown in Figure 3. This is given by: 

logioNec = logioNe logioNe < 5.5 

logioN^c = logioNe x (1 + {log^N^ - 5.5) x 0.04) log^oNe > 5.5 



Where N^c is the corrected electron size, is the original size and E is the 
shower energy in PeV as described above. This increases the electron size by 
around ~ 7% near energies of lOPeV. This simple parameterization provides 
a much better fit to the shower curves in Figure 3 as shown by the filled 
symbols identified as 'saturation corrected data'. We use this corrected electron 
size through the rest of this paper. Importantly this change is not large, but 
illustrates the power of correlating apparently redundant measurements. These 
can be used to explore systematic problems with various measured quantities. 



4 Mass Determination 

In this analysis we adopt a slightly different philosophy from previoTis work 
by using combinations of measured parameters to derive an incident parti- 
cle mass, ^4, on a shower by shower basis. With the information presented 
above we form two estimates for the mass, one from the location of X^ax 
and the fitted shower energy and another from the muon and electron sizes in 
combination with the fitted energy. 

4-1 Mass Functions 

Although mass functions of arbitrary complexity can in principle be derived 
from analysis of numerical simulations, these are not always transparent as to 
the underlying physical processes which produce the sensitivity to mass from 
a particular combination of parameters. Here we adopt the approach of using 
a simple physical model to produce a mass function which is tested against 
simulations to determine the accuracy of this function over some range of 
energies. We first make some simple assumptions about the variation of shower 
parameters for singly charged particles with energy: 

Xmax = -'^O + X logio{E) 
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Where Xo, Xm, -/V/iO, are constants and E is the energy in PeV. Typical 
values for the scaling constants are Xm = SOg/crn^ and a = 0.87. We assume 
simple superposition to define how these will vary for heavier particles of mass 
A: 

Xmax = Xo + XmX logw{E / A) 

ijLg^logioiA) + ax logw{E/A) 



Where jig — logio{Ni^/Ni^o) and the first of these equations is the same as A. 2. 
We can rearrange these to give an estimate for the mass, A: 



Where logio{A)x is the mass from the measurement of Xmax and logiQ{A)g is 
the mass from the measurement of muon size at ground level. At fixed energy, 
deeper showers (larger Xmax) give a smaller mass and muon rich showers 
(larger jig) give a larger mass. Another feature of these mass estimate equations 
is that they have an opposite sign dependence on the shower energy. Hence 
errors in energy assignment produce opposing errors in the mass estimates. 
In principle these estimates can be used to make independent mass estimates 
from the same showers. Equation 1 results from the location of the hadronic 
shower maximum and equation 2 arises from the energy of the pions decaying 
into muons in the shower. At the same total energy, proton showers have 
typically higher energy pions than iron showers so heavier nuclei produce more 
muons from pion decays in the shower development. 

There is however another factor which effects the muon size at ground level. 
This is slightly dependent on X^ax in a similar but much weaker fashion to the 
electron size dependency discussed above and shown in Figure 3. Fluctuations 
in shower development produce larger muon sizes at ground level when Xmax 
is larger for the same initial energy. Although this is not a large effect, the 
mass estimate described in equation 2 is sensitive to small fluctuations in /ig 
because of the small value of (1 — a) in the denominator. These fluctuations 
are strongly correlated with the electron size at ground level and their effect 
can be significantly reduced by introducing a correction term based on Cg = 
ln{Nec/Neo), the normalized electron size at ground level. The overall mass 
function derived from muon measurements using this correction is given by: 



logw{A)^ = {Xq - Xmax)/Xm + logw{E) 
logiQ{A)g ^{/ig-ax logio{E))/ (1 - a) 



(2) 
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logio{A)g ^{jig-ax logio{E) - S{eg))/{1 - a) 



(3) 



Where the correction S{eg) is apphed to reduce fluctuations in muon size 
arising from variations in the shower development. This correction is discussed 
in detail in Appendix C. 

4-2 Correlations of Mass Estimates 

To investigate the consistency of these mass estimates we can look at the 
variations of one type of mass estimate versus the other. This is shown in 
Figure 4 for some of the energy bins studied. The correlation between these 
estimates over a range of energies provides confidence that the estimates are 
roughly correct. This correlation has some small contribution from the com- 
mon energy term in equations 1 and 2, this is not large since a fluctuation 
in energy of 30% only produces a change in either mass function of ~ 0.1. 
The dominant contribution to this correlation is from the rapid variation of 
electron size measured at the ground with X^ax depth as shown in Figure 3. 
Although clearly //.^ plays a role, since there are more muons from showers 
produced by heavier primaries, the fluctuations of showers at some flxed pri- 
mary energy and mass provide wide enough variations in X^ax to make the 
electron size effect predominate. In some sense the correlations in Figure 4 are 
a straightforward result of the relationships shown in Figure 3. 

4-3 Comparison with Simulations 

To further test these distributions and methods for estimating mass we can 
compare these results with the expectations from simulations. These calcu- 
lations must take account of shower properties, detector properties, and the 
analysis cuts used. A simulation program which uses the results of air showers 
produced by the CORSIKA 4.5 code as an input for the shower characteristics. 
The simulation effort was concentrated on the DICE measurements of X^^ax: 
all the details of the fluctuations in the quantities directly affecting these were 
included. The electron sizes and muon sizes at ground level are included but 
the fluctuations in these arc parameterized in a simplified fashion to speed 
up the calculations. These simulated showers impact on the detector site at 
random locations and angles near the vertical in the simulation. The output 
of the simulation is a data file which can be subjected to exactly the same 
cuts and analysis procedures applied to the actual data. For the analysis dis- 
cussed here the input primary mass composition used in the simulation is a 
constant 30% H, 30% He, 10%CNO, 10% NeMgSi, and 20% Fe. As discussed 
later, there is no evidence for signiflcant variation in this composition across 
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Fig. 4. Cross correlation plots of the mass determinations from X^nax and the muon 
size for events in energy bins similar to Figure 2. The size of each box is linearly 
proportional to the number of events, concentration of the events are expected along 
the dashed lines if the correlation is good. The mass functions are discussed in the 
text. 

the range of energies in this work. The simulation is also generated with a 
spectral shape which approximates the variation in slope through the knee 
region, and includes details of the contributions from shower geometry, focal 
plane pixel size, etc. The simulation can also be used to determine the overall 
effective geometrical factor of the instrument. This was used in combination 
with the instrument livetime to calculate the normalization of the data in 
Figure 1. 

We can compare the shape of the X^ax mass estimates from DICE with the 
predictions from the simulation. This provides a test of the overall consistency 
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Fig. 5. (Left) Xjnax mass distribution RMS widths compared with simulations across 
all measured energy bins. The right hand axis is labelled with the equivalent RMS 
width in X^ax- (Right) Detailed Xmax mass distribution at 4PeV compared with 
simulations discussed in the text. The top axis is labelled with the corresponding 
values of X„iax using equation A. 2 

of the method used and can exclude certain extreme composition possibilities. 
The left hand panel of Figure 5 shows the RMS value of logio{A) compared 
with expectations from simulations for mixed composition discussed above, 
hght (p+a), and heavy (A>4) composition. The comparison of the complete 
distribution for the Xmax mass at 4PeV is shown in the right panel of Figure 
5 as data points with errors. This also shows the expected distributions for a 
mixed composition (line), light (dashes), and heavy (dots). Although this is 
not the method of choice for the most accurate composition determination, 
these data alone seem to exclude a pure Fe composition over all of this energy 
range. Unlike other methods which use the distribution of Cherenkov light at 
ground level, the mass determination of DICE docs not produce a long tail for 
logio{A)x in the negative region. This is because DICE intrinsically rejects all 
showers with an apparent Xmax deeper than ~ 750g/cm^ by geometry. This 
provides some truncation of the distribution produced of the mass function 
which is incorporated in the simulation. 

With distributions such as shown in the right hand panel of Figure 5 the mass 
resolution of the DICE logio{A)x measurements can be estimated. These are 
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Fig. 6. The RMS resolution derived from simulations for single events for logio{A) 
derived from the measurement of Xj^ax by DICE. The dashed curve is for proton 
events, the solid line for iron nuclei. 



shown in Figure 6 for incident protons and iron nuclei as a function of parti- 
cle energy. The ordinate shows the RMS resolution of the DICE experiment 
for the logio{A)x of a single event. The improvement in the resolution with 
increasing energy arises from a combination of larger light yield and an im- 
proving Xmax resolution with the better event geometry in DICE of deeper 
showers. The mass resolution of proton showers at energies above a few PeV 
is close to the intrinsic limit for a mass based on X^ax determination. Even 
if the depth of shower maximum and the shower energy is known perfectly, 
the mass resolution in logiQ(A) for protons is limited to a ~0.8 because of 
intrinsic fluctuations in the shower development. The mass resolution for iron 
showers is significantly better than for protons because the intrinsic fluctua- 
tions of showers from these nuclei are significantly smaller. It is sobering to 
realize that iron is separated by only ~ 2(7 from an incoming proton signal for 
single events. 
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5 Composition Determination 

Using these methods a mass estimate for each event can be derived from (i) 
^max -I (ii) muon and electron size, and (iii) a simple linear combination of 
the two, called logiQ{A)c. This last method has the advantage of somewhat in- 
creased resolution because the methods are to some extent independent. Each 
event is assigned values for logiQ{A)xi logio{A)g, and logio{A)c and collected 
in bins of energy. The results for the mean logio{A) values versus energy are 
shown as the upper panel in Figure 7. The mass error is based on statistics 
and the energy error is the estimate of the systematic error in the energy scale 
of the fits discussed earlier. All of the mass estimates show variations with 
energy across the knee region. These data are tabulated in Appendix B. 

Another method for determining composition has been used (see e.g. Wat- 
son [27]) which attempts to separate a sample of almost pure proton+a events 
by applying a simple cut to the mass distribution, such as shown in Figure 
5, which excludes heavy particles. The simulation can then be used to make 
corrections to the samples produced by estimating the amount of 'leakage' of 
heavier events into the (proton+a) sample. The results of an analysis on these 
data which uses a cut to extract the events where logiQ{A) < are shown in 
the lower panel of Figure 7 for the DICE X^nax results. Similar analyses have 
not been made for the ground size distributions because simplifications in 
the fiuctuations for these data make this method unreliable. The errors on 
the (p+a) fraction are based on simple counting statistics, the uncertainty is 
larger for this analysis method because fewer events pass the data cuts. There 
is a tendency of these results to a higher fraction of (p+a) at the higher en- 
ergies although they are also consistent with constant composition across this 
region. The values of the (p+a) fraction are tabulated in Appendix B. 



6 Discussion 

By using several parameters measured for each shower a mass estimate for each 
incoming primary particle can be extracted which shows remarkable consis- 
tency between values derived from Xmax and from the /j/e method. The direct 
observation of the expected shape and absolute normalization of showers in 
this energy range (Figure 3) is a crucial test for the absence of systematic 
errors in the methods of mass estimation. The widths of the distributions 
shown in Figure 5 are consistent with the results shown Figure 7. These re- 
sults essentially exclude an all-iron composition at 10 PeV, this conclusion 
was reached by an earlier analysis of DICE. The results suggest a tendency to 
lighter composition at high energies, although those from the ratio method are 
consistent with constant composition. A composite mean logio{A)c using both 
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Fig. 7. (Top) Results for mean logiQ{A) measurements from this work versus energy, 
the model curves are discussed in the text. (Bottom) Derived fraction of (p+a)/all 
for the Xmax mass determination. 
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estimates shows a decreasing mass across the energy range of the knee region. 
For comparison, two models of possible cosmic ray mass variation near the 
knee are plotted in Figure 7. The dotted curve is the model given by Swordy 
1995, where a source spectrum oc is introduced above the knee region and 
the energy dependent escape from the Galaxy is assumed to reach a plateau 
near the knee region. This model provides an increase in mass across the knee 
of a size which is larger than for models based simply on particle rigidity. 
The lower dashed line, with which the data seem more consistent, is a model 
which has similar low energy behavior as the previous model but which in- 
troduces a proton source which compensates for the lost flux above the SNR 
acceleration cutoff rigidity, assumed here to be lO-'^^V. This is similar to the 
suggestion of Protheroe and Szabo 1992 [28] where an extragalactic compo- 
nent of particles produced by AGNs provide the cosmic rays at high energy. 
These data seem more consistent with this latter model, where cosmic rays 
become progressively lighter across the knee region. However, they do not ex- 
clude the possibility that the composition is more or less constant across this 
energy range. A sudden change in composition to becoming predominantly 
iron nuclei {logio{A) = 1.75) seems strongly excluded. The variations in the 
energy spectrum of Figure 1 near 8 PeV and the irregularity in the behavior 
of logio{A) rms' near 4PeV shown in Figure 5 lend some support to the sin- 
gle source model of the knee region discussed by Erlykin and Wolfendale[29]. 
However, Figure 7 does not show the expected increase in < logio{A) > of 
~0.4 across the knee region predicted by this model. 

Contrary to conventional wisdom these results do not support a simple 'rigidity 

steepening' which would lead to a steady increase in mass across the knee 
region. If the cosmic ray abundances below the knee are provided with a simple 
steepening from a spectrum of E~'^''^^ to and a composition below the 

knee is similar to that described by Swordy[22], the size of the effect should 
be an increase of ~0.2 in the mean value of logio{A) across the knee. The 
combination of experimental measurements discussed here have the sensitivity 
to detect an increase of this order if it were present. 

The mass resolution of these methods is poor when compared with direct 
techniques. The la (RMS) widths of the resolution in logio{A) for nuclei have 
typical size a ~ 0.8. Hence protons which fluctuate upward by ~ 2a can 
be misidentified as iron nuclei. In this measurement environment any meth- 
ods which produce a lower intrinsic systematic error in the technique are 
extremely valuable. The philosophy behind DICE is to try to reduce these 
types of error by decreasing the reliance on Monte Carlo simulations of air- 
showers in the interpretation of the data. The method of directly measuring 
Xmax from geometry avoids steps associated with inferences about the loca- 
tion of Xmax from other variables such as the density distribution of particles 
at ground level which are interpreted with simulations. Only the probability 
of Cherenkov emission as a function of angle with the shower axis is needed to 
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estimate X^ax with DICE. This is relatively insensitive to the details of the 
simulation. 



The use of Xmax as a variable for composition determination is desirable for 
other reasons. The mass estimate derived from this variable is less sensitive 
to uncertainties in the characterization of shower parameters. To explore the 
sensitivity of the results to the analysis procedures a simple test can be made 
on the data set by exploring the trend of the derived masses for different values 
of the analysis parameters. The results of this test are shown in Figure 8. Here 
the upper panel shows two lines which are the trend of the data, as shown in 
Figure 7, under two assumptions for the value of the elongation rate, X^n in 
equation A. 2. The solid line uses X^^ — SO g/cm^ in the analysis, the dashed 
line uses X^ — 75g/cm^. The smaller elongation rate increases the apparent 
lightening of composition by an amount of ~0.1 in < LogiQ{A) >, giving 
an estimate of the systematic uncertainties in the analysis of the DICE Xmax 
data. The lower panel shows a similar test of the n/e mass derived from ground 
particle densities. Here the solid line represents the superposition model where 
a — 0.87 in equation 3, the dashed line is an analysis where a — 0.82. This 
mass estimate is clearly sensitive to the precise value of a used since this 
change makes < Logio{A) > apparently increase by an amount ~0.2 across 
this energy region. The systematic error in the fi/e mass is therefore large 
since simulations typically give values of a in the range ~ 0.8 — 0.9. The above 
analyses illustrate a crucial property of multiparameter measurements in this 
type of work: the requirement for consistency from a variety of apparently 
redundant measurements can directly assess the validity of the result. In the 
case discussed here the mean mass estimates from the different techniques 
must agree since these are the same events, consequently an analysis where 
Xjn — Ibg/cw? and a — 0.82 is excluded by this consistency check. 



The instrumentation of ground based cosmic-ray air shower detectors have im- 
proved significantly over the past decade or so. The statistical quahty of these 
measurements has improved to the point where the intrinsic fluctuations in 
the shower developments are becoming the limiting factor in the determina- 
tion of primary composition. The new challenge emerging in the derivation of 
composition from measurements of airshowers is the identification and quan- 
tification of systematic errors in aspects of the measurements themselves or in 
the simulations used for analysis. By making use of multiparameter measure- 
ments of the same showers we can directly explore our level of understanding 
of these issues and provide a reliable determination of cosmic-ray composition 
in this energy region. The results presented here are an attempt to take some 
initial steps in this direction. 
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Fig. 8. (Top) The sensitivity of the Xmax mass measurement to the elongation rate 
assumed. The solid line shows the trend of the data with the rate 80g/cm? as used 
in this analysis, the dashed line shows the trend of an analysis using Ybg/cm^. 
(Bottom) The sensitivity of the /i/e mass measurement to the analysis parameters. 
The solid line shows the trend of the data analyzed with a superposition model 
where N^^ oc £''^•87 as used here. The dashed line shows the trend produced if the 
analysis assumes oc E^'^"^ . 
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Appendices 



A Appendix of Fitting Functions 

A.l Cherenkov Light Distribution Function 

This function which is used to fit the Cherenkov hght sizes at the DICE 
detectors is derived from simulations using CORSIKA 4.5 with VENUS [19] 
and is vahd for shower core distances larger than 75m from each DICE cluster. 
It depends on particle energy, mass and event geometry as follows: 

Ch = CoE{PeVyA-'e-^'-^"'^ (A.l) 



where 



Co = 1-89 X lO'' 
7 = 1.144 + 0.0905 X logw{A) 
/3 = 0.0161 X (1 - 0.128 X logw{A) 

+0.124 X logio{E{PeV)) x (1 + 0.322 x logio{A))) 
6 = 0.186 



Here Ch is the number of Cherenkov photons collected in the effective area 
(3.75m2) of a DICE detector. The total energy of the shower in PeV is E{PeV), 
A is the mass of the incident particle, and r(m) is the perpendicular distance 
to the shower core from the DICE detector. Note the amount of light varies 
almost linearly with energy (7) since the Cherenkov size is close to a calorime- 
teric measurement of the shower energy because nearly all of the particles pro- 
duced in the shower give Cherenkov hght. The light intensity also decreases 
with increasing mass for the same total energy (e). The light intensity falls 
with distance from the shower core, for low mass particles near IPeV the light 
intensity falls with a scale distance of ~62m {1//3). For iron nuclei at similar 
energies the scale distance is ~80m. 
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A. 2 Xjnax Fitting Function 

The following equation is the fitting function used for the primary mass and 
energy dependence for the location of X^ax' 

Xmax = + X {logioiE) - logio{A)) (A.2) 

where 

Xo = 560g/cm^ 
Xm = SO^/cm' 



Where X^ax is the location of shower maximum, E and A are the incident 
particle energy in PeV and particle mass (in nuclear mass units) respectively. 

A. 3 Gaisser Shower Size Function 

The function used for the average shower development is adapted from Gaisser 
[26] to calculate the expected electron size at ground level, Xq, for a shower 
of total energy E with a maximum at a depth Xjnax'- 

Ne(Xo) = X E(PeV) X l)e^(Xo/(X^„, - A))(^+^)e-^°/^ (A.3) 

Here Ag — 6.02 x 10^ is a constant, E{PeV) is the shower energy in PeV, 

Xq = 860g/cm'^ is the average detector atmospheric depth, p — X^ax/^ — 1 
and A = 70g/cm'^ is the hadronic interaction length. This assumes simple 
scaling in the Gaisser function, i.e. 5*0 = 0.045 and e = 0.074 GeV. The lower 
right panel of Figure 3 shows the effect of the scaling violation suggested by 
Gaisser[26] as a dashed curve. 
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B Appendix of Tabulated Results 

This Appendix contains tables of the results of this investigation. 
Table B.l 

Table of X^ngx versus energy from DICE 



Energy Bin 


Median Energy 




(PeV) 


(PeV) 


gjcvr? 


0.2-0.4 


0.26 ± 0.05 


470.7 ± 0.6 


0.4-0.7 


0.50 ± 0.10 


476.3 ± 0.8 


0.7-1.5 


0.91 ± 0.18 


497.7 ± 1.0 


1.5-3.0 


1.92 ± 0.38 


529.6 ± 1.6 


3.0-5.0 


3.66 ± 0.73 


562.8 ± 3.3 


5.0-8.0 


6.03 ± 1.21 


587.4 ± 5.5 


8.0-15.0 


10.09 ± 2.02 


610.3 ± 10.0 



Table B.2 

Table of mass estimates from DICE, ground sizes from CAS A/MIA and from a 
combination of the two methods. Also given is the result of the (p-|-a)/all fraction 



estimate, discussed in the text. 


Energy 


Combined 




n/e 


{p+a)/all 


(PeV) 


< /0510(A) > 


< logio{A) > 


< logiQ{A) > 




0.2 


0.61 ± 0.01 


0.64 ± 0.01 


0.76 ± 0.01 


0.609 ± 0.011 


0.5 


0.71 ± 0.01 


0.82 ± 0.01 


0.73 ± 0.02 


0.560 ± 0.012 


1.0 


0.70 ± 0.01 


0.79 ± 0.01 


0.68 ± 0.02 


0.531 ± 0.014 


2.0 


0.63 ± 0.02 


0.70 ± 0.02 


0.58 ± 0.03 


0.494 ± 0.022 


4.0 


0.53 ± 0.04 


0.54 ± 0.04 


0.53 ± 0.05 


0.577 ± 0.045 


6.0 


0.47 ± 0.06 


0.51 ± 0.07 


0.45 ± 0.09 


0.633 ± 0.082 


10.0 


0.34 ± 0.10 


0.40 ± 0.12 


0.39 ± 0.12 


0.674 ± 0.112 
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Table B.3 



The All Particle Fluxes of Cosmic Rays from DICE 


lotai x^arxicie rL/iiergy (^ijev ) 


r LUX [TLl ST S VjiC V ) 


r lux H/rror (^it^^m sr s vjCV ) ) 


U.ooox lU 


U.ODi X lU 


n Kr\ 1 n~12 
U.oUX iU 


U.44 ( X lU 


n one \/ 1 n~io 

U.ZUoX iU 


n in \/ 1 n~i2 
U.oUX iU 


U.ODZ X lU 


n 1 1 t^s/ 1 n— 10 

U. 1 iO X lU 


n 1 a V 1 n— 12 

U. io X lU 


U. (UoX iU 


n f^i Ox/ 1 n~ll 
U.Diy X lU 


n 1 1 \/ 1 n~12 
U.il X iU 


U.oyi X iU 


n Qoo \y 1 n^ll 


U.DO X lU 


n 1 1 \/ 1 n7 


n 1 K7\/ 1 n^ll 

U.ID ( X lU 


u.oy X lu 


A 1 /1 1 \/ 1 n7 

U.141 X lU 


U.yUoX iU 


n OA \y 1 n~i3 

U.Z4X lU 


U.i /oX iU 


U.4yDX iU 


n 1 cc V/ 1 n~13 
U.iOX iU 


n 00/1 1 n7 


U.ZDo X iU 


n ofis/ 1 n— 14 
u.yox lu 


U.ZoZ X lU 


n 1 /10s/ 1 n— 12 

U.14ZX iU 


n f?i s/ 1 n— 14 

U.Di X lU 


U.OOO X lU 


n 7^/Lv 1 n— 13 

U. ( D4X lU 




0.447x10^ 


0.333x10^13 


0.23x10-1^ 


0.562x10^ 


0.163x10^13 


0.14x10-1"^ 


0.708x10^ 


0.711x10^1'' 


0.82x10-15 


0.891x10^ 


0.557x10-1^ 


0.64x10-15 


0.112x10^ 


0.207x10-1^ 


0.35x10-15 


0.141 xlO^ 


0.999x10-15 


0.21x10-15 


0.178x10^ 


0.430x10-15 


0.12x10-15 


0.224x10^ 


0.113x10-15 


0.57x10-16 


0.282x10^ 


0.112x10-15 


0.50x10-16 
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C Appendix detailing the muon mass function 



As has been mentioned in the main text the mass derived from the muon size, 
equation 3, is susceptible to variations introduced by fluctuations from the 

location of shower maximum. Although these are not as large as associated 
with the electron size shown in Figure 3, they produce poor resolution in the 
muon mass function unless corrected. A quantitative measure of this effect is 
shown in Figure C.l 

The normalized muon mass size {logiQ{Nfj_/ E{PeV)) is derived from the MIA 
muon size and the DICE Cherenkov energy, E(PeV). This size is clearly cor- 
related with Xmax depth at energies near 1 PeV as shown. This is a weaker 
version of the effect shown for electrons in Figure 3. It presumably arises from 
more muons either decaying or reaching the end of the range for the showers 
with maxima higher from the ground. Although at a fixed X^ax the muon 
mass function will perform as desired, the presence of fluctuations smears the 
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Fig. C.l. The variation of the normahzed loglO(muon size/energy) with depth of 
shower maximum (Xmax) around 1 PeV is shown as the filled symbols. There is a 

clear correlation between the depth of shower maximum and the normalized muon 
size. The lower dotted line is a simulation of this effect for masses A < 4, the upper 
dashed line is a simulation for yl > 4 using the same composition discussed in the 
text. 
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possible mass resolution because of this effect. Indeed as the relationship be- 
tween the simulation curves (dashed lines) on Figure C.l and the measured 
data show, in general as the shower penetrates more deeply the muon con- 
tent appears more 'proton-hke' as expected. To produce a more accurate mass 
function a variable must be introduced to offset the effect of shower depth 
shown in the Figure. Since electron size is closely correlated with X^ax, as 
shown in Figure 3, this can be used to make a correction and provide a mass 
function from the ground sizes which is completely independent from the Xjn^x 
measurement of DICE and therefore provide a separate mass estimate for each 
shower. 

Figure C.l shows the basis for this correction. For showers near 1 PeV the av- 
erage muon size, normalized by the Cherenkov energy, is plotted against the 
normalized electron size. The variations in these sizes come predominantly 
from the fluctuations in X^ax between showers. Figure C.l represents on av- 
erage how the muon size changes with electron size due to these fluctuations. 
The dashed line is a linear flt to these data with a slope of 0.5. The electron 
size correction to the muon size mass measure is therefore of the form: 

S{eg) = 0.5 X {eg - 1.14 x logioE{PeV)) (C.l) 

Where Cg — logio{Nec/ Ngo) as before, and where the second term corrects 
for the average electron size increase with energy. This effectively provides a 
correction for the muon size growth with X^ax shown in Figure C.l. 
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Fig. C.2. The variation of the normahzed loglO(muon size/energy) with normaUzed 
loglO(electroii size/energy) near 1 PeV. The hne is a fit with a slope of 0.5, this is 
discussed further in the text 
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